function [x1_cup,y_cu,diagSim,y_cup,Zfixed] = solutionExtendedPath(Z0,x_cu,y_tN,N,setupEPer)

setupEPer.dispOn = 0;
f0               = DSGEforesight_N(reshape(Z0,setupEPer.ny+setupEPer.mx,N),x_cu,y_tN,N,setupEPer); 

if setupEPer.fixedPointSolver == 1
    [Zfixed,fval,diff,iter,errorMes] = FixedPointSolver(Z0,x_cu,y_tN,N,setupEPer);
    usedSolver = 1;
elseif setupEPer.fixedPointSolver == 2
    [Zfixed,fval,diff,iter,errorMes] = FixedPointSolver_deltaOpt(Z0,x_cu,y_tN,N,setupEPer);
    usedSolver = 2;
elseif setupEPer.fixedPointSolver == 3
    % The unweighted NLS optimization
    MaxIter  = 500*2;
    weightOn = 0;
    [Zfixed,fval,diff,iter,Q,errorMes] = FixedPointSolver_LM(Z0,x_cu,y_tN,N,MaxIter,weightOn,setupEPer);
    usedSolver = 3;
elseif setupEPer.fixedPointSolver == 4
    % The weighted NLS optimization
    MaxIter  = 500;
    weightOn = 1;
    [Zfixed,fval,diff,iter,Q,errorMes] = FixedPointSolver_LM(Z0,x_cu,y_tN,N,MaxIter,weightOn,setupEPer);
    usedSolver = 4;
end
if errorMes ~= 0 && setupEPer.fixedPointSolver <= 3
    if setupEPer.dispOn == 1
        disp('Baseline fixed-point solver failed. New attempt by a mix of the routines')
    end
    % Improve starting values for delta_opt using a few iterations by the
    % LM minimizer
    MaxIterLM = 10;
    if fval'*fval < f0'*f0
        ZstartLM = Zfixed;
    else
        ZstartLM = Z0;
    end
    weightOn = 0;
    Z0 = FixedPointSolver_LM(ZstartLM,x_cu,y_tN,N,MaxIterLM,weightOn,setupEPer);
    
    % Try with delta_opt
    MaxIterDefault = setupEPer.MaxIter;
    setupEPer.MaxIter = 200;
    [Zfixed,fval,diff,iter,errorMesBackup2] = FixedPointSolver_deltaOpt(Z0,x_cu,y_tN,N,setupEPer);
    setupEPer.MaxIter = MaxIterDefault;
    usedSolver = 2;
    if errorMesBackup2 == 1
        Z0 = Zfixed;
        MaxIterLM = 1000;
        weightOn  = 0;
        [Zfixed,fval,diff,iter,Q,errorMesBackup3] = FixedPointSolver_LM(Z0,x_cu,y_tN,N,MaxIterLM,weightOn,setupEPer);
        usedSolver = 3;
    end
end
Zfixed                    = real(Zfixed);
y_cu(:,1)                 = Zfixed(1:setupEPer.ny,1);
y_cup(:,1)                = Zfixed(1:setupEPer.ny,2);
x1_cup(1:setupEPer.nx1,1) = [Zfixed(setupEPer.ny+1:setupEPer.ny+setupEPer.mx)';
                            y_cu(1:setupEPer.myx,1)'];
diagSim                   = [max(abs(fval)) max(abs(fval(1:setupEPer.ny+setupEPer.mx,1))) iter usedSolver];


end

